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ABSTRACT 


Patterns of broad-scale plant species richness are thought to be largely determined by (1) variation in energy and water 
availability among sampling units (species energy hypothesis), (2) habitat and topographic heterogeneity within sampling units 
(spatial heterogeneity hypothesis). and (3) regional differences in geographic configuration and history (regional effects 
hypothesis). However, lack of taxonomic and distribution data, particularly for tropical regions, has impeded assessments of 
the relative importance of these three hypotheses. We used a large botanical database to estimate the pattern of relative 
vascular plant richness across the western Neotropics and regression models to measure the extent to which this estimated 
pattern supported predictions from each of the above three hypotheses. Variation in plant richness across three major 
paleophysiographic regions (northwest South America, southern Central America, and northern Central America) was 
primarily predicted by the spatial heterogeneity hypothesis, with secondary contributions from the species energy hypothesis 
and, to a lesser extent, the regional effects hypothesis. Regression models that incorporated the relative contributions of all 
three hypotheses predicted peaks of relative species richness mostly in topographically complex areas (e.g., Sierra Madre de 
Chiapas, Cordillera de Tilarán, Cordillera de Talamanca, Panama’s Cordillera Central, the Andes, and the Venezuelan 
Guayana); relatively low richness in central Mexico and Yucatán, Los Llanos of Venezuela, and in the Gran Chaco region of 
Bolivia, Paraguay, and Argentina; and a richness trough in lowland Amazonia relative to southern Central America, the Andes, 
and the Venezuelan Guayana. We discussed the contrast between our results and previous assessments that found plant 
richness to be primarily determined by the species energy hypothesis and predicted different patterns of plant richness across 
the western Neotropics. 

Key words: Broad-scale species richness, Neotropics, regional effects, spatial heterogeneity, species energy, vascular 
plants. 


Patterns of spatial variation in the number of 
species co-occurring within broad geographic areas 
were discussed by naturalists of the 18th and 19th 
centuries (von Humboldt, 1808; Darwin, 1862; 
Wallace, 1878; see Hawkins, 2001), and more 
rigorous attempts to quantify these patterns began 
during the last century (e.g.. Wulff, 1935; Simpson, 
1964; see Mutke & Barthlott, 2005) along with the 
formulation of numerous hypotheses to explain them 
(Hutchinson, 1959; Pianka, 1966; Huston, 1979, 
1994; Rohde, 1992; Palmer, 1994; Rosenzweig, 
1995; Schemske, 2002; Willig et al., 2003). Yet, 
understanding what determines broad-scale spatial 
patterns of species richness still remains a central 
issue in ecology and biogeography (Gaston, 2000; 
Ricklefs, 2004; Pennisi, 2005). Recent work in this 


area has focused on a few prominent hypotheses, 


notably the species energy (SE), spatial heterogeneity 
(SH), and regional effects (RE) hypotheses. 

The SE, broadly defined, holds that gradients of 
energy and water availability across sampling units 
create and maintain species-richness patterns (Hutch- 
inson, 1959; Pianka, 1966; Brown, 1981; O’Brien, 
1995) through a variety of 
nonexclusive mechanisms (Hawkins et al., 2003; 
Currie et al., 2004; Evans et al., 2005; Clarke & 
Gaston, 2006) including the effects of these gradients 
on extinction. (Wright, 1983), evolutionary rates 
(Rohde. 1992; Evans & Gaston, 2005). and species 
filtering according to physiological requirements (von 
Humboldt, 1808; Turner et al., 1987). The SH posits 


that variation in elevation and habitat within sampling 


1993; Rosenzweig, 


units increases richness, again through various 


nonexclusive mechanisms, including species sorting 
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among habitats (Rosenzweig, 1995) or isolated 
patches of similar habitats (Simpson, 1964) and the 
effect of spatial heterogeneity on diversification rates 
(Simpson, 1964; Jetz et al., 2004). The RE proposes 
that differences in richness between sampling units 
with similar environments but located in different 
regions (i.e., diversity anomalies) result from regional 
disparities in history and geographic configuration 
that, in turn, cause differences in the occurrence of 
particular lineages, time available for diversification, 
and rates of species production and extinction 
(Pianka, 1966; Latham & Ricklefs, 1993; Schluter & 
Ricklefs, 1993; Qian & Ricklefs, 2000). 

Many ecologists seem to agree that some combina- 
tion of the three hypotheses above provides the best, 
currently available, explanation of broad-scale spatial 
patterns of species richness (Rahbek & Graves, 2001; 
Ricklefs et al., 2004; Field et al., 2005; Kreft & Jetz, 
2007). However, progress in establishing the relative 
importance of these hypotheses as explanations of 
broad-scale patterns of plant species richness has 
been hindered by the scarcity of taxonomically and 
spatially resolved data sets (Whittaker et al., 2005), 
particularly for tropical regions (Frodin, 2001; but see 
Kier et al., 2005). Thus, most studies are focused on a 
few relatively well-known life forms (e.g., Currie & 
Paquin, 1987; O'Brien, 1998; Field et al., 2005) or 
taxonomic groups (e.g., Bjorholm et al., 2005). Other 
studies are based on floras, checklists, and other 
literature sources that have little information on 
sampling effort (e.g., Kreft & Jetz, 2007). Yet, failure 
to account for sampling effort can alter our perception 
of spatial patterns of plant richness (Nelson et al., 
1990; Parnell et al., 2003). Here, we expanded on 
previous efforts to overcome these problems (Jiménez 
et al., 2009) by assessing the relative importance of 
the SE, SH. and RE in determining broad-scale 
patterns of vascular plant species richness across the 
western Neotropies. In addition, we used results on 
the relative importance of different hypotheses to 
predict relative plant species richness in areas of the 
western Neotropics with poorly sampled floras. 


METHODS 


To describe the spatial pattern of vascular plant 
species richness across the western Neotropics, we 
used 755,401 georeferenced herbarium specimen 
records from the Tropicos database (<http://www. 
tropicos.org/>, Missouri Botanical Garden), repre- 
senting 48,264 plant species collected in central 
Mexico, Belize, Guatemala, Honduras, Nicaragua, 
Costa Rica, Panama, French Guiana, Guyana, Vene- 
zuela, Colombia, Ecuador, Peru, Bolivia, Uruguay, 
and northern Argentina. We mapped these data on a 


Behrmann cylindrical equal-area projection of the 
study area and used rarefaction (Gotelli & Colwell, 
2001) to estimate relative species richness in 100 X 
100 km sampling units as the expected number of 
species in a random sample of size n from a set of N 
herbarium specimen records. To choose the value of 
n at which relative richness was measured using 
rarefaction, we balanced the benefits of increasing n 
in terms of increased precision of relative richness 
estimates against the corresponding costs in terms of 
decreased sample size (Fig. 1; Jiménez et al., 2009). 
We judged n — 500 herbarium specimen records to be 
a reasonable compromise and retained 255 sampling 
units with at least 500 specimen records for further 
analysis (Fig. 2A, B). 

Admittedly, rarefaction is unlikely to account for all 
error in estimates of relative richness due to 
differences in sampling effort, because the set of N 
herbarium specimens collected in a single sampling 
unit is unlikely to be a random sample from the 
individuals of all vascular plants occurring in that 
sampling unit (Jiménez et al., 2009). Therefore, our 
estimates of relative plant richness contain potentially 
large measurement errors that, nonetheless, can be 
reasonably subsumed in the error term of statistical 
models representing different hypotheses. Such mea- 
surement errors may affect the conclusions from our 
analysis mainly in proportion to its correlation with 
the explanatory variables relevant to the three 
hypotheses of interest (Jiménez et al., 2009). As a 
starting point, we assumed such correlation was 
negligible. 

We used a set of regression models to simulta- 
neously measure the extent to which the estimated 
pattern of relative plant richness across the western 
Neotropics (Fig. 2A) supported predictions from the 
SE. SH, and RE. For each regression model, we 
derived predictions from one or more hypotheses 
about the sign of regression coefficients relating 
explanatory variables to relative plant richness. Some 
regression models included higher-order terms (i.e., 
interaction or quadratic terms) that may account for 
small portions of the variation in the response variable 
when the range of the relevant explanatory variable is 
limited. Therefore, when the coefficients of higher- 
order terms were not statistically significant, we 
examined the performance of reduced models with 
no higher-order terms in an attempt to trade-off a 
tolerable bias for increased precision (Chatterjee & 
Hadi, 1988). In these reduced models, the predicted 
sign for regression coefficients may depend on the 
range and central tendency of explanatory variables. 
We derived predictions for the sign of regression 
coefficients in reduced models based on the range of 
explanatory variables in our sample (Fig. 3) and the 
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Trade-off between precision and sample size. —A. Pearson’s correlation coefficients (line and circles) and their 


95% confidence intervals (dotted lines) for the relationship between richness estimated by rarefaction as the expected number 
of species in 5000 specimen records and richness estimated by rarefaction at various other numbers of specimen records (in 
the abscissa): 4000, 3000, 2000, 1000, 500, and 200. The sample to estimate all these correlations was a set of 41 sampling 
units of 100 X 100 km, each with at least 5000 herbarium specimen records. Also shown is the potential sample size (line and 
triangles): the number of sampling units in our study area that would be available for analysis if richness estimates based on 
rarefaction at a given number of specimen records were deemed acceptable. By deciding to estimate richness at 500 specimen 
records we adopted a somewhat conservative approach, trading off increases in potential sample size above 255 sampling units 
for relatively precise richness estimates. —B. Relationship between relative plant richness estimated by rarefaction at 5000 
herbarium specimen records (in the abscissa) and 500 herbarium specimen records (in the ordinate) across 41 sampling units 


of 100 X 100 km, each with at least 5000 specimen records. 


values of regression coefficients obtained by previous 
studies (O’Brien, 1998; Francis & Currie, 2003; Field 
et al., 2005). 

The SE was represented by three models thought to 
explain a major portion of broad-scale spatial variance 
in plant richness worldwide (O’Brien, 1998; Francis & 
Currie, 2003; Field et al., 2005; Table 1). Importantly, 
we selected data sources to reflect the original 
formulation of each model as closely as possible. Thus, 
for SE models 1 and 2 (Francis & Currie, 2003; 
Table 1) we obtained annual potential evapotranspira- 
tion and water deficit from Ahn and Tateishi (1994), 
while for SE model 3 (O'Brien, 1998; Field et al., 2005; 
Table 1) we calculated minimum monthly potential 
evapotranspiration using mean monthly temperature 
and Thornthwaite's formula (Thornthwaite, 1948; see 
details in Jiménez et al., 2009). We obtained data for 
mean annual temperature, mean monthly temperature, 
and annual precipitation from WorldClim (<http:// 
www.worldclim.org/>; Hijmans et al., 2005). 

The SH was represented by three regression models 
based on hypothesized effects on plant richness of 
within-sampling-unit spatial variation in elevation 
(Currie & Paquin, 1987; O'Brien et al., 2000; Kreft & 
Jetz, 2007), climate (Currie & Paquin, 1987; Linder, 
2003). and soil (Linder, 2003; Tuomisto et al., 2003; 


Gentry, 1988; Table 1). We calculated variance in 
elevation within each sampling unit from elevation 
data at a resolution of 90 X 90m (Shuttle Radar 
USGS, 2004). We 


spatial variation in climate within sampling units 


Topography Mission, derived 


from 19 variables measuring various aspects of 
temperature and precipitation at a resolution of 30 
arcseconds, obtained from WorldClim. Principal 
component analysis on these 19 variables captured 
89.4% of the total variation in the first three principal 
components. Because variance within sampling units 
in these three principal components was highly 
correlated (Pearson's r > 0.75; P « 0.05), we used 
only the variance within sampling units in the first 
principal component as an explanatory variable. We 
caleulated variance in available water capacity and 
soil carbon density within each sampling unit from 
data at a 5 X 5 min. resolution (Global Soil Data 
Task). 

We also considered regression models that com- 
bined terms representing the SE and SH. One model 
was the Interim General Model second-generation 
(IGM2) (O'Brien et al., 2000; Field et al, 2005; 
Table 1), regarded as one of the best working models 
to explain broad-scale patterns of woody plant 
richness (O'Brien et al., 2000; Field et al., 2005; 


A B 


WU 435 - 461 12827 - 39247 
EE 428-434 8325 - 12826 
EH 425-427 5725 - 8324 
EN 421-424 3902-5724 
Ml 417-420 3403-3901 

412-416 2726-3402 

| 406-411 2081-2725 

|402-405 1823-2080 

|396-401  1614- 1822 
380-395 1407-1613 
EE 384-389 1155 - 1406 
EX 378-383 1049-1154 
EN 370-377 961-1048 
EN 362 - 369 842 - 960 
EN 355 - 361 785 - 841 
MMM 347 - 354 719 - 784 
EN 339 - 346 658 - 718 
EN 327 - 338 595 - 657 
WME 308 - 326 545 - 594 
EN 201 - 307 500 - 544 


20°N 
C D 


WE 418-450 963-2120 
BS 412-417 852 - 962 
ES 405 - 411 782 - 851 
ES 400 - 404 731 - 781 
396 - 399 690 - 730 
| 392 - 395 639 - 689 
| 388 - 391 595 - 638 
| 385 - 387 536 - 594 
| 383 - 384 480 - 535 
380 - 382 419 - 479 
EH 377 - 379 357 - 418 
ES 373 - 376 287 - 356 
EH 369 - 372 197 - 286 
ES 365 - 368 121 - 196 
EN 361 - 364 48 - 120 
EN 356 - 360 -59 - 47 
E 349 - 355 -151 - -60 
EN 339-348 -223 - -152 
"WE 327-338 -300 - -224 
EH 234-326 -553 - -301 


10°N 


10°S) 


20°S 


l [ | i [] [! 
TOW 400W sowo Bowo TOW Sowo 50W 

Figure 2. Relative plant richness and botanical sampling effort in 255 sampling units of 100 X 100 km distributed across the western Neotropics. Each color represents an interval of five 
percentiles, from high (red) to low (blue). —A. Observed richness measured as the average number of species found in 1000 rarefied samples of 500 herbarium specimen records. —B. Number of 
herbarium specimen records. —C. Predicted relative richness based on the fit of the data in A to the nonspatial version of Interim General Model second-generation (IGM2) with regional effects 
(Table 2). Prediction is restricted to the sample range of explanatory variables. —D. Woody plant species richness as predicted by the global coefficients for the IGM2 (Field et al., 2005). No 
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Figure 3. Distribution of explanatory variables across the western Neotropics (gray bars). and our sample of the study area 
(bars with black outline). At the top of each histogram. a line and a dot represent the range and mean, respectively. for the 
western Neotropics (gray), our sample of the study area (black), and samples used to generate species energy (SE) models and 
the globally specified Interim General Model second-generation (IGM2) (dashed line). The summary statistics for values of 
mean annual temperature (A), water deficit (D), and potential evapotranspiration (B) in the samples used to generate SE 
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Clarke & Gaston, 2006). While the IGM2 was 
specifically developed for woody plants, its concep- 
tual basis is general and allows application to all 
plants, at least as an approximate representation of the 
SE and SH. To replicate as closely as possible the 
original formulation of the IGM2, we estimated 
elevation range for each sampling unit using the 
GTOPO30 arcsecond elevation data set (USGS, 2004) 
resampled to 10 X 10 km. We constructed additional 
models combining terms representing the SE and 
SH that yielded significant (P < 0.05) regression 
coefficients that were consistent with the respective 
predictions in all previous regression models. We 
examined all combinations of these explanatory 
variables, except when they were highly correlated. 
For brevity, we presented only the models with more 
support Akaike 
Information Criterion corrected (AICc) for small 


empirical as measured by the 


sample sizes, recommended when the ratio of 
sampling units to model parameters is less than 40 
(Burnham & Anderson, 2002). Given the procedure 
used to construct these latter models, we refer to them 
as ad hoc models (Table 1). 

The RE was represented by regression models that 
also included terms representing both the SE and SH 
(Table 2). This approach follows from the idea that 
regional effects account for variation in species 
richness that remains after the effect of environmental 
conditions of the sampling units has been accounted 
for (Schluter & Ricklefs, 1993). We added regional 
effects to regression models using dummy variables 
(Draper & Smith, 1998) to code intercepts and slopes 
for each of three major paleophysiographic regions of 
(Graham, 1997): 


through. northern. Costa Rica (N Central. America), 


the Neotropics central Mexico 
southern Costa Rica and Panama (S Central America), 
and northern South America (N South America). There 
are at least two well-known major differences in the 
history of these paleophysiographic regions (Graham, 
1997; Burnham & Graham, 1999): (1) N South 
America was isolated from N Central America for 
lens of millions of years, and (2) unlike N South 
America and N Central America, S Central America 
emerged from the sea just a few million years ago. In 
addition, the three regions were probably differential- 
ly influenced by orogenic activity, and by late 
Cenozoic fluctuations in climate and sea level 
associated with glacial advances and retreats that 
intensified during the Quaternary (Graham, 1997; 


Bumham & Graham, 1999). Following previous 


approaches to examine regional effects (Schulter & 
Ricklefs, 1993; Ricklefs et al., 2004), we predicted 
differences among regions but not the direction of 
such differences. Variables representing the RE were 
removed from regression models if they did not reduce 
AlCe values, following model simplification proce- 
dures suggested by Crawley (2002). 

To confront the estimated pattern of relative plant 
richness across the western Neotropics (Fig. 2A) 
against the models representing different hypotheses, 
we used quantile regression through the median (Cade 
et al., 2005), implemented using R package quantreg 
(Koenker, 2005, 2008; R Development Core Team, 
2006). We used regression through the median rather 
than least squares, because the response variable 
measures plant richness in an ordinal scale and 
meaningful hypotheses about variables in an ordinal 
scale focus on order statistics such as the median 
(Wolman, 2006). We gauged the extent to which the 
data supported different regression models examining 
the statistical significance of regression coefficients 
and their concordance with the predictions derived 
from the respective hypotheses. We also assessed 
support for different regression models using AlCc. 
We calculated model fit in terms of R’, the proportion 
of the sum of absolute deviations from the median of 
the response variable that is accounted for by a 
regression model (Cade et al., 2005). 

All regression models incorporated a covariate 
(logarithm of area) to account for reduced area in 
sampling units intersecting the border of the spatial 
extent of the study. Pairwise correlations among 
(Table 3) that no predictors 
included together in a single regression model were 


predictors showed 
highly correlated and, thus, collinearity was unlikely 
to be an issue. Furthermore, collinearity due to 
higher-order regression terms was alleviated by mean- 
centering the explanatory variables of models with 
higher-order terms (Quinn & Keough, 2003). We 
tested for spatial autocorrelation in regression resid- 
uals using a permutation test for Moran's I (Fortin & 
Dale, 2005), implemented with R package spdep 
(Bivand, 2008). When spatial dependence in the 
residuals was significant, we used spatial eigenvector 
mapping (SEVM; Dormann et al., 2007) to construct 
spatial regression models that estimated the relation- 
ship between relative species richness and explana- 
tory variables in the absence of spatial autocorrelation 
in regression residuals. We used a forward-selection 
procedure to first include in quantile regression 


= 


models come from Francis and Currie (2003). Summary statistics for values of annual mean precipitation (C) and minimum 
monthly potential evapotranspiration (E) in the samples used by Field et al. (2005) to generate the globally specified IGM2 
come from O’Brien (1998). and those for elevation range (F) are from O'Brien et al. (2000). 
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Table 1. Results from median regression models representing the species energy (SE) and spatial heterogeneity 


(SH) hypotheses. 


Regression coefficient, Regression coefficient, 
Models and variables Pred nonspatial spatial 
Model 1 SE, R? — 0.039 (0.212) 
Log10 area 55.34 (23.89) 67.15 (20.26) 
Water deficit 0.02 (0.02) 0.02 (0.02) 
Mean annual temperature + 0.13 (0.11) 0.24 (0.07)* 
Water deficit X mean annual temperature 2X 103*15:2€.10:3) ] 2€ 10/*432X 105 $ 
Eigenvector 3 119.08 (40.63) 
Eigenvector 2 210.64 (51.58) 
Eigenvector 5 109.42 (43.95) 
Eigenvector 4 138.61 (40.18) 
Eigenvector 12 186.20 (55.10) 
Eigenvector 6 164.99 (52.71) 
Model 1 SE reduced, R? — 0.037 (0.221) 
Log10 area 51.200 (24.81) 44.58 (23.11) 
Water deficit 0.021 (0.02) 0.02 (0.01) 
Mean annual temperature + 0.151 (0.12) 0.12 (0.07) 
Eigenvector 2 176.65 (53.93)* 
Eigenvector 3 131.52 (40.62) 
Eigenvector 5 98.36 (41.72)* 
Eigenvector 12 193.03 (52.09) 
Eigenvector 6 164.89 (49.37) 
Eigenvector 15 158.76 (47.53)° 
Model 2 SE, R! — 0.084 (0.177) 
Log10 area 47.80 (24.79 42.72 (24.06) 
Water deficit 0.01 (0.02) 0.002 (0.01) 
Annual potential evapotranspiration + 0.05 (0.02)? 0.05 (0.02)* 
Annual potential evapotranspiration? 3 X 10 *(8 X 10 4+ 3X 10 *(3 X 10 *» 
Eigenvector 2 194.63 (49.93) 
Eigenvector 12 219.29 (52.03)° 
Eigenvector | 7.78 (51.64) 
Eigenvector 3 93.58 (44.48) 
Model 3 SE, R? — 0.109 (0.200) 
Logl0 area 63.52 (12.68) 42.15 (17.37) 
Annual precipitation + 0.01 (2 X 10 ?y 0.01 (0.003) 
Min. monthly potential evapotranspiration + 0.15 (0.13) 0.04. (0.15) 
Min. monthly potential evapotranspiration? 0.01 (3 X 10 ?p 0.02 (0.01)* 
Eigenvector 2 164.31 (47.71) 
Eigenvector 8 117.62 (51.81) 
Eigenvector 12 234.95 (57.95) 
Model 1 SH, R! — 0.143 (0.298) 
Log10 area 42.59 (15.44 37.97 (10.86) 
Log10 variance in elevation + 16.477 (2.4)° 19.84 (2.11)* 
Eigenvector 2 158.58 (29.357 
Eigenvector 3 139.58 (27.87) 
Eigenvector 12 125.16 (45.66) 
Eigenvector 11 100.16 (33.61) 
Eigenvector 10 137.96 (33.53) 
Model 2 SH, R! — 0.132 (0.253) 
Log10 area 44.80 (24.91) 33.88 (22.81) 
Log10 variance in climate PC 1 + 66.718 (6.37) 60.75 (9.49)° 
Eigenvector 2 173.27 (46.80) 
Eigenvector 3 124.51 (41.92 
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Table 1. 


Continued. 


Regression coefficient, 


Regression coefficient, 


Models and variables Pred nonspatial spatial 
Eigenvector 12 147.93 (40.46) 
Eigenvector 15 101.07 (41.73) 

Model 3 SH, R? — 0.051 (0.220) 
Log10 area 74.34 (23.89) 53.16 (17.08) 
Log10 variance in available water capacity + 10.21 (5.52) 12.78 (3.89) 
Log10 variance in soil carbon density + 20.43 (8.33)? 22.46 (6.00) 
Eigenvector 2 126.86 (44.30) 
Eigenvector 3 102.54 (37.80) 
Eigenvector 7 121.81 (37.62) 
Eigenvector 12 218.07 (48.86) 
Eigenvector 6 146.83 (40.59) 
IGM2, R! — 0.218 (0.274) 
LoglO area 50.79 (8.43) 55.5] (9.93)° 
Annual precipitation + 0.011 (2 X 10 °% 0.01 (0.003)* 
Min. monthly potential evapotranspiration + 0.298 (0.11)? 0.28 (0.12)* 
Min. monthly potential evapotranspiration? 0.003 (3 X 10 ?) 0.01 (0.003) 
Ln elevation range + 15.037 (1.99)° 15.06 (2.01) 
Eigenvector 2 117.08 (41.17) 
Eigenvector 8 90.23 (42.97) 
Eigenvector 3 93.55 (39.41) 
Model 1 ad hoc, R? — 0.254 (0.298) 
LoglO area 52.70 (8.83) 59.35 (8.40)* 
Annual precipitation + 0.01 (0.002)° 0.01 (0.003)* 
Min. monthly potential evapotranspiration + 0.29 (0.10) 0.27 (0.12) 
Min. monthly potential evapotranspiration? 0.01 (0.003)* 0.002 (0.003) 
LoglO variance in elevation + 18.77 (1.987 19.64 (2.49)° 
Eigenvector 2 110.38 (40.88)? 
Eigenvector 3 109.42 (38.23) 
Eigenvector 8 109.61 (43.72)? 
Model 2 ad hoc, R — 0.258 (0.300) 
Log10 area 49.93 (6.29) 64.49 (7.37)° 
Annual precipitation * 0.01 (7 X 10 *)* 0.01 (0.003) 
Min. monthly potential evapotranspiration * 0.280 (0.10)? 0.25 (0.11)* 
Min. monthly potential evapotranspiration? 0.01 (0.002)? 0.003 (0.003) 
Log10 variance in soil carbon density + 5.19 (3.38) 4.54 (3.59) 
Log10 variance in elevation + 18.477 (1.9) 18.65 (2.35)° 
Eigenvector 2 94.64. (38.38 
Eigenvector 3 98.03 (34.85) 
Eigenvector 8 102.67 (40.44) 
Model 3 ad hoc, R? — 0.254 (0.300) 
Log10 area 43.43 (10.82) 60.95 (16.24) 
Annual precipitation * 0.01 (0.003) 0.01 (0.002 
Potential evapotranspiration + 0.001 (0.02) 0.03 (0.02) 
Potential evapotranspiration? 2x 10 *(-2 X 10 ‘F PO bx 0" 4) 
Log10 variance in soil carbon density + 8.289 (3.90)* 6.65 (4.45) 
Log10 variance in elevation + 16.025 (1.62)° 16.88 (2.07)° 
Eigenvector 2 136.09 (39.45) 
Eigenvector 3 122.15 (31.55y 
Eigenvector 8 90.70 (40.97 
Model 4 ad hoc, R? — (0.303) 
Log10 area 47.47 (15.19) 
Potential evapotranspiration + 0.02 (0.01) 
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Table 1. 


Continued. 


Regression coefficient, Regression coefficient, 


Models and variables Pred nonspatial spatial 
Potential evapotranspiration? 2:26: TOA KTO 
Log10 variance in elevation + 16.19 (2.047 


Eigenvector 2 


174.42 (34.07)° 


The response variable was relative plant species richness. The first column shows names of models and explanatory 
variables, as well as goodness of fit (R?) for nonspatial regression models first and for spatial models in parentheses. The 
column labeled “Pred” shows the predicted signs of regression coefficients according to the respective hypothesis. The next 
two columns show regression coefficients for nonspatial and spatial models with standard errors in parentheses and statistical 
significance coded as: *, P < 0.05: °”, P < 0.01: °, P < 0.001. PC, principal component; IGM2, Interim General Model second- 


generation; Log10, logarithm base 10; Ln, natural logarithm. 


models those spatial eigenvectors that most reduced 
the spatial autocorrelation in the regression residuals 
(see Jiménez et al., 2009). Given that the interpreta- 
tion of differences between coefficients derived from 
spatial and nonspatial regression models remains 
controversial (Dorman et al., 2007), we provided 
results from both types of models. 

To predict relative plant richness across the 
western Neotropics, including areas that have been 
only sparingly collected, we used the predicted values 
from the regression models best supported by the 
data—that is, models that yielded statistically signif- 
icant regression coefficients, were consistent with the 
respective predictions, and had the lowest AICe 
values. To avoid extrapolation, we predicted richness 
only in sampling units that fell within the range of 
explanatory variables in our sample of the western 
Neotropics (Fig. 3). 


RESULTS 
NONSPATIAL MODELS 


The performance of nonspatial models representing 
a single hypothesis indicated that the SH had more 
empirical support than the SE. All three models 
representing the SH yielded significant regression 
coefficients consistent with the respective predictions, 
while only two out of three models representing the SE 
did so (Table 1). The AICe values for SH models 1 
and 2 were notably lower than those for other models 
representing a single hypothesis, differing by > 17 
(Fig. 4) and indicating large differences in empirical 
support. SH models | and 2 also accounted for a 
higher proportion of the variation in relative species 
richness (R! — 0.143 and 0.132, respectively) than 
other models (Table 1). Thus, the relative perfor- 
mance of SH models 1 and 2 suggested a major role 
for variance (within sampling units) in elevation and 
determinants of plant 


climate, respectively, as 


richness across the study area. 


The performance of models combining terms 
representing the SE and SH was, with no exception, 
superior to that of models representing a single 
hypothesis. Each of these four models yielded 
significant regression coefficients consistent with 
predictions from both the SE and SH (Table 1). AICe 
values for the IGM2 and ad hoc models 1 to 3 were 
substantially smaller than those for models represent- 
ing an individual hypothesis, differing by > 41 
(Fig. 4). In addition, the IGM2 and ad hoc models 1 to 
3 accounted for more variation in relative richness (R! 
> 0.21) than models representing a single hypothesis 
(Table 1). The higher performance of ad hoc models 1 
to 3 compared to models representing a single 
hypothesis and to the IGM2 could be partly due to 
over-fitting, given the ad hoc procedure used to build 
the former models (see Methods). However, we used 
the IGM2 as an a priori model and, therefore, at least 
in this case, the higher performance of a model that 
combines hypotheses relative to models that represent 
a single hypothesis indicates the complementary 
nature of the SE and SH. 

All models that included terms representing the RE 
yielded significant regression coefficients consistent 
with predictions from both the SE and the SH. Each of 
these models also yielded significant coefficients 
representing regional effects (Table 2), but we had no 
predictions about their sign (see Methods). Models 
representing the RE did best overall, with lower AICc 
values differing by > 7 from models combining terms 
representing the SE and SH with no regional effects 
(Fig. 4), although respective differences in model fit 
were slight (Tables 1, 2). In general, the improvement 
in model performance accomplished by adding regional 
effects was smaller than that achieved by combining 
terms representing the SE and SH. This finding 
suggests that the RE, as represented in this study, 
does not complement the SE and SH to the same degree 
that the latter two hypotheses complement each other. 

Nonetheless, adding regional effects to models 
representing the SE and SH decreased AlCe values 
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substantially. For example, adding regional effects to 
the IGM2 decreased the AlCe value by 23. The 
resulting model suggested that, after controlling for 
other environmental variables, relative richness was 
higher in S Central America than in N South America, 
and higher in N South America than in N Central 
America, across the range of minimum monthly 
potential evapotranspiration and annual precipitation 
in our sample (Fig. 5A, B). In addition, maximum 
plant richness was attained at higher values along the 
axis of minimum monthly potential evapotranspiration 
in N Central America than in S Central America and 
N South America (Fig. 5A). Relative plant richness 
increased faster with elevation range (within sampling 
units) in N Central America than in S Central America 
and N South America. Relative richness in sampling 
units with narrow elevation ranges was lower in N 
Central America than in S Central America and N 
South America, but the difference in relative richness 
decreased as elevation range increased (Fig. 5C). 
There were differences among ad hoc models in the 
strength and significance of regional effects, but all 
models revealed that N Central America, compared 
with S Central America and N South America, had a 
lower intercept and a steeper slope relating relative 
richness to variation in elevation within sampling 


units (Table 2). 


SPATIAL MODELS 


Similar to the results from nonspatial models, the 
performance of spatial models representing a single 
hypothesis indicated more empirical support for the SH 
than the SE. All three models representing the SH 
yielded significant regression coefficients consistent 
with the respective predictions (Table 1). However. SH 
model 3 also yielded a significant negative regression 
coefficient for variance in potential available water 
capacity, contrary to the respective prediction. Two of 
three models representing the SE, SE models 2 and 3, 
yielded significant regression coefficients consistent 
with the respective predictions, while SE model 1 
yielded a significant negative regression coefficient for 
mean annual temperature, contrary to the respective 
prediction (Table 1). Among spatial models represent- 
ing a single hypothesis, SH models 1 and 2 explained 
more variation in the response variable (R? — 0.298 
and 0.253, respectively; Table 1) and had substantially 
lower AICe values than the rest (Fig. 4), corroborating 
results from nonspatial models and suggesting a 
primary role for variance in elevation or climate, 
within sampling units, as determinants of plant 
richness across the study region. 

Each of the five spatial models combining terms 
representing the SE and SH yielded statistically 


significant regression coefficients that were consistent 
with the predictions derived from both hypotheses. In 
no case was there a statistically significant coefficient 
opposite to any prediction (Table 1). These models 
performed better than those representing an individ- 
ual hypothesis, with a notable exception: the spatial 
version of SH model 1, which described a positive 
relationship between variance in elevation. within 
sampling units and relative plant richness, explained 
more variation in the response variable and had a 
substantially lower AICe value than the IGM2, in both 
metrics similar to ad hoc models | to 4 (Table 1; 
Fig. 4). This exception strengthens the previous 
suggestion that variance in elevation within sampling 
units, a variable representing the SH, was a primary 
determinant of plant richness across the study region. 
It also suggests that spatial eigenvectors accounted for 
variation in relative plant richness that correlated with 
variables representing the SE, but not for variation in 
relative plant richness that correlated with variables 
representing the SH. 

All spatial models including terms representing the 
RE yielded statistically significant regression coeffi- 
cients consistent with predictions from both the SE 
and SH, as well as significant coefficients represent- 
ing regional effects. In only one case was there a 
statistically significant coefficient opposite to a 
prediction: ad hoe model 3 with regional effects 
yielded a negative coefficient for potential evapo- 
transpiration (Table 2). Spatial models including 
regional effects explained only a slightly higher 
proportion of the variation in the response variable 
than models including terms representing the SE and 
SH only (Tables 1, 2), and their AICe values were not 
consistently lower than those of other models (Fig. 4). 
This result contrasts with the respective comparison 
for nonspatial models and suggests that spatial 
eigenvectors accounted for variation in relative plant 
richness that correlated with variables representing 
the RE. Nonetheless, a spatial model including 
regional effects yielded the lowest AlCe (ad hoc 
model 2; Fig. 4) and the highest proportion of 
explained variation in relative species richness (R' 
— 0.326; Table 2). This model revealed similar 
regional effects to those described by nonspatial 
models, whereby N Central America compared to 
other paleophysiographic regions had a lower inter- 
cept and a steeper slope relating relative richness to 
variance in elevation within sampling units (Table 2; 


Fig. 5D-F). 


PREDICTED PLANT RICHNESS MAPS 


The models best supported by the data yielded 
similar patterns of predicted relative species richness 
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Table 2. Results from median regression models combining the species energy, spatial heterogeneity, and regional 
effects hypotheses. 


Regression coefficient, Regression coefficient, 
Models and variables Pred nonspatial spatial 
IGM2, R? — 0.263 (0.297) 
LoglO area 55.54 (13.49) 64.749 (18.32) 
Annual precipitation + 0.01 (0.002)° 0.01 (0.003)° 
Min. monthly potential evapotranspiration + 0.15 (0.14) 0.07 (0.13) 
Min. monthly potential evapotranspiration? 0.01 (0.004)? 0.01 (0.004)? 
Ln elevation range 10.39 (1.95) 10.55 (1.91) 
CA 112.35 (21.65) 111.42 (30.07)* 
SCA 23.63 (11.96) 29.22 (9.24) 
Min. monthly potential evapotranspiration X NCA 0.87 (0.28 0.70 (0.35)* 
Ln elevation range X NCA 16.00 (3.23) 15.15 (4.34) 
Eigenvector 3 140.92 (31.76) 
Eigenvector 2 14.53 (53.74) 
Model 1 ad hoc, R? — 0.284 (0.308) 
LoglO area 50.87 (7.61) 65.41 (7.34) 
Annual precipitation + 0.01 (0.002)° 0.01 (0.002)° 
Min. monthly potential evapotranspiration + 0.14 (0.13 0.001 (0.12) 
Min. monthly potential evapotranspiration? 0.01 (0.004)° 0.01 (0.003) 
Log10 variance in elevation + 14.73 (2.49% 14.13 (2.18) 
NCA 70.49 (22.79) 60.06 (28.04)* 
SCA 59.05 (50.46) 102.91 (53.73) 
Annual precipitation X SCA 0.03 (0.01)° 0.04 (0.01)? 
Min. monthly potential evapotranspiration XX NCA 0.72 (0.25) 0.69 (0.25) 
Log10 variance in elevation X NCA 14.06 (4.67)? 10.95 (5.61)* 
Log10 variance in elevation X SCA 24.58 (10.52) 33.43 (11.63)? 
Eigenvector 3 153.32 (29.98) 
Model 2 ad hoc, R? — 0.277 (0.326) 
LoglO area 48.65 (9.44) 50.08 (13.57) 
Annual precipitation + 0.01 (0.002)° 0.01 (0.003) 
Min. monthly potential evapotranspiration * 0.31 (0.11) 0.33 (0.12) 
Min. monthly potential evapotranspiration? 0.01 (0.003) 0.01 (0.003) 
Log10 variance in soil carbon density + 10.15 (3.72)? 3.77 (3.40) 
Log10 variance in elevation + 15.36 (2.42)° 15.36 (2.49) 
NCA 111.96 (30.75) 86.68 (36.12) 
Annual precipitation X NCA 0.01 (0.01)* 2 X 10 * (0.01) 
Log10 variance in elevation X NCA 17.84 (5.40) 17.90 {6.11} 
Eigenvector 8 98.63 (35.58) 
Eigenvector 2 149.03 (40.37) 
Eigenvector 10 112.11 (31.40y 
Model 3 ad hoc, R? — 0.291 (0.304) 
Log10 area 48.25 (8.73) 44.42 (11.41) 
Annual precipitation * 0.01 (0.003)* 0.01 (0.003)* 
Annual potential evapotranspiration * 0.04. (0.02) 0.06 (0.02) 
Annual potential evapotranspiration? 38 x3107* (5€ 1975) 3210. * (d 2€ 10-5y 
Log10 variance in soil carbon density + 10.81 (3.85)? 13.95 (4.55)? 
Log10 variance in elevation + 12.22 (1.95) 11.24 (2.45y 
NCA 56.45 (17.41) 62.95 (22.46) 
SCA 25.56 (9.52)? 26.78 (12.07 
Annual potential evapotranspiration XX NCA 0.14 (0.05)° 0.17 (0.04) 
Log10 variance in elevation X NCA 10.74 (3.59) 12.21 (4.55)? 
Eigenvector 7 97.57 (42.37)* 
Model 4 ad hoc, R? — (0.272) 
Log10 area 51.32 (13.73)* 
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Table 2. Continued. 
Regression coefficient, Regression coefficient, 
Models and variables Pred nonspatial spatial 

Annual potential evapotranspiration + 0.02 (0.02) 
Annual potential evapotranspiration? 2 X 10 *(2 X 10 * 
Log10 variance in elevation + 15.11 (1.94)° 
NCA 17.78 (7.04) 
SCA 32.58 (7.357 
Annual potential evapotranspiration X NCA 0.07 (0.05) 
Eigenvector 3 184.55 (35.16) 
Eigenvector 2 83.90 (33.53) 


The response variable was relative plant species richness. The first column shows names of models and explanatory variables, 
as well as goodness of fit (R?) for nonspatial regression models first and for spatial models in parentheses. The column labeled 
“Pred” shows the predicted signs of regression coefficients according to the respective hypothesis. The next two columns show 
regression coefficients for nonspatial and spatial models with standard errors in parentheses and statistical significance coded 
as: *, P < 0.05; ^, P < 0.01; *, P < 0.001. IGM2, Interim General Model second-generation; Log10, logarithm base 10; Ln, 
natural logarithm: NCA, dummy variable for N Central America; SCA, dummy variable for S Central America. 


across the western Neotropics (Figs. 2C, 6A-D). 
Generally, the highest species richness was predicted 
in topographically complex areas such as the 
mountains to the east of the Isthmus of Tehuantepec 
including the Sierra Madre de Chiapas, the mountain 
ranges extending from the Cordillera de Tilarán 
southeast along the Cordillera de Talamanca into 
Panama's Cordillera Central, the Andes, and the 
Venezuelan Guayana. Areas predicted to have highest 
richness formed a longitudinally broad band in 
Colombia, encompassing the Chocó region and all 
three Andean cordilleras, and included both Andean 
cordilleras in Ecuador. These areas were largely 
restricted to the eastern Andes in Peru and Bolivia 
according to most models, with the exception of ad hoc 
model 4, which predicted high plant richness in the 


Table 3. 
and spatial heterogeneity hypotheses. 


western portion of the Andes in Peru (Fig. 6C), 
presumably because this latter model did not include 
any variable explicitly measuring water availability. 
All models predicted relatively low richness in central 
Mexico and Yucatán, Los Llanos of Venezuela, and in 
the Gran Chaco region of Bolivia, Paraguay, and 
Argentina. All models also predicted a species 
richness trough in lowland Amazonia relative to S 
Central America, the Andes, and the Venezuelan 
Guayana. These predicted patterns of relative plant 
richness should be seen in the light of important 
differences between observed and predicted richness 
values. Specifically, even the models best supported 
by the data accounted for relatively small portions of 
the variation in the response variable (R! — 0.258- 
0.326; Tables 1, 2). 


Pearson's correlation coefficients among variables used in regression models representing the species energy 


Area PET WD Temp Precip mPET Elev Range Soil Pawc PCAI 
PET 0.352 1.000 
WD 0.234 0.121 1.000 
Temp 0.301 0.489 0.155 1.000 
Precip 0.184 0.448 0.518 0.392 1.000 
mPET 0.321 0.611 0.123 0.775 0.536 1.000 
Elev 0.137 0.070 0.282 0.596 0.121 0.382 1.000 
Range 0.172 0.061 0.261 0.563 0.118 0.375 0.951 1.000 
Soil 0.256 0.084 0.050 0.055 0.145 0.121 0.069 0.061 1.000 
Pawe 0.095 0.058 0.129 0.083 0.013 0.060 | 0.055 0.033 0.545 1.000 
PCAI 0.073 0.020 0.220 0.444. 0.031 0.20601 0.838 0.799 0.161 0.112 1.000 
Rich 0.166 0.149 0.146 0.094. 0.301 0.105 0.399 0.346 0.130 0.008 — 0.401 


Area, log10 area; PET, potential evapotranspiration: WD, water deficit; Temp. annual mean temperature; Precip, annual 
precipitation; mPET, minimum monthly potential evapotranspiration: Elev, logl0 variance in elevation: Range. In elevation 
range; Soil, log10 variance in soil carbon density; Pawc, log10 variance in available water capacity: PCA1, log10 variance in 
climate’s first principal component; Rich, relative plant richness measured by rarefaction at 500 herbarium specimen records. 
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Figure 4. Support for models representing the species energy (SE), spatial heterogeneity (SH), and regional effects (RE) 
hypotheses estimated by differences in Akaike Information Criterion corrected (AICc) for small sample sizes (delta AICc) 
between each model and the model performing best (1.e.. the model with the lowest AICc). The horizontal dotted lines separate 
models representing a single hypothesis (top), models combining terms representing the SE and SH (middle), and models 
combining terms representing the three hypotheses (bottom). Triangles represent spatial models and circles nonspatial models. 
Filled symbols represent regression models that yielded statistically significant regression coefficients consistent with the 
respective predictions. Open symbols represent regression models that yielded no statistically significant regression 
coefficients consistent with the respective predictions. Crosses mark regression models that yielded statistically significant 
regression coefficients that were inconsistent with the respective predictions. When higher-order terms (interaction or 
quadratic terms) were not statistically significant, we examined the performance of reduced models with no higher-order terms. 
Therefore, some models are represented by two circles or triangles. 


150 


Disci ssion plant richness to explanatory variables were supported 
at least once (potential evapotranspiration and its 
Our results supported most a priori predictions square, minimum monthly potential evapotranspiration 


based on previous studies about determinants of plant and its square, annual precipitation, within-sampling- 


richness (Currie & Paquin, 1987; Gentry, 1988; 
O'Brien, 1998; O'Brien et al., 2000; Francis & Currie, 
2003; Linder, 2003; Tuomisto et al., 2003; Field et al.. 
2005; Kreft & Jetz, 2007). Nine of the 13 predictions 


regarding the sign of regression coefficients relating 


unit range in elevation, within-sampling-unit variance 
in elevation, within-sampling-unit variance in climate, 
and within-sampling-unit variance in soil carbon 
density); three predictions were opposed once (annual 
mean temperature, potential evapotranspiration, and 
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Figure 5. Plots of partial residuals (Draper & Smith, 1998) for models that incorporate the species energy (SE). spatial 
heterogeneity (SH), and regional effects (RE) hypotheses: —A-—C. Nonspatial version of Interim General Model second- 
generation (IGM2) with regional effects (Table 2). —D-—F. Spatial version of ad hoc model 2 with regional effects (see Table 2). 
These plots show the estimated relationship between relative plant richness and the variable in the abscissa when all other 
variables in the model are statistically controlled except for regional effects. Open circles (and dashed lines) represent N 
Central America, open triangles (and solid lines) represent S Central America, and filled circles (and solid lines) represent N 
South America. 
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within-sampling-unit variance in available water 
capacity); and two predictions were neither supported 
nor opposed (water deficit and the interaction between 
water deficit and annual mean temperature) by 
statistically significant results (Tables 1, 2). We also 
found statistically significant differences in relative 
plant richness among similar environments located in 
different paleophysiographic regions (Table 2; Fig. 5) 
defined according to previous work (Graham, 1997; 
Burnham & Graham, 1999). These findings suggest that 
we used a reasonable set of hypotheses and corre- 
sponding models to examine the relative importance of 
the major determinants of broad-scale plant richness 
across the western Neotropics. 

The performance of nonspatial and spatial regression 
models, measured by the consistency of regression 
coefficients with a priori predictions and by AICe and 
R' values, indicated that explanatory variables repre- 
senting the SH were the primary determinants of plant 
richness across the western Neotropics, with comple- 
mentary contributions from variables representing the 
SE and, to a lesser extent, the RE. In particular, 
variance in elevation and variance in climate within 
sampling units, representing the SH, were the main 
predictors of the estimated pattern of relative plant 
richness across the study area. The models performing 
best overall included variance in elevation and not 
variance in climate, but distinguishing which of these 
two variables was more important was difficult because 
they were highly correlated (Table 3). A third variable 
representing the SH, variance in soil carbon density 
within sampling units, was included in some of the best 
performing models, but its role was secondary to that of 
variance in elevation or climate. Annual precipitation, 
minimum monthly potential evapotranspiration, and 
potential evapotranspiration, representing the SE, most 
effectively improved the performance of models 
representing the SH. Here, again, determining whether 
minimum monthly potential evapotranspiration or 
potential evapotranspiration was more important was 
difficult because they were correlated (Table 3). 
Finally, variables representing the RE sometimes 
increased the performance of models combining terms 
representing the SH and SE. Most frequent among 
these were variables indicating a lower intercept for N 
Central America than for other regions, a higher 
intercept for S Central America than for other regions, 
and higher slopes relating plant richness to spatial 
heterogeneity within sampling units (variance and 
range in elevation) and to water (annual precipitation) 
and energy availability (potential evapotranspiration 
and minimum monthly potential evapotranspiration) in 
N Central America than in other regions. 

One of the main findings emerging from our 
analysis was more empirical support for the SH than 


for the SE, consistent with a similar study of broad- 
scale plant richness across northwestern South 
America that found at least as much support for the 
SH as for the SE (Jiménez et al., 2009). This result 
would seem at odds with previous work indicating that 
variables representing the SE are more important 
determinants of broad-scale patterns of plant richness 
than those representing the SH (Currie & Paquin, 
1987; O'Brien et al., 2000; Hawkins et al., 2003; 
Bjorholm et al., 2005; Field et al., 2005; Moser et al., 
2005; Kreft et al., 2006; Kreft & Jetz, 2007). We are 
aware of only one earlier plant study (Pausas et al., 
2003) documenting a primary role for the SH outside 
our study region. Below, we explore potential 
explanations for the differences between our results 
and those from previous studies, related to the 
characteristics of the response variable and the 
distribution of the explanatory variables. 

Geographie patterns of plant richness measured in 
small sampling units (e.g. = lha. plots) are 
sometimes considered broad scale (e.g.. Hawkins et 
al.. 2003), but the importance of different hypotheses 
can be contingent on sampling unit size. For example, 
the relationship between plant richness and variables 
representing the SE may be most evident when 
measured across large sampling units, while other 
factors may exhibit greater heterogeneity and thus be 
more important determinants of spatial variation in 
richness at smaller grains (Whittaker et al., 2001). 
This would seem to be at least as much of an issue for 
variables representing the SH. Specifically, the SH 
proposes that spatial heterogeneity fosters species 
coexistence across habitats and isolated patches of 
similar habitat, or that it accelerates speciation rates 
by increasing opportunities for isolation and ecolog- 
ical divergence. Both of these effects are likely 
increasingly opposed by dispersal as sampling unit 
size decreases (cf. Moser et al., 2005). Thus, the 
relative importance of the SH may increase with 
sampling unit size, as suggested by studies of bird 
richness (Rahbek & Graves, 2001; van Rensburg et 
al., 2002; Hulbert & Haskell, 2003). 

Regarding procedures to estimate the response 
variable, previous studies estimated richness by 
overlapping geographic range maps (e.g., O’Brien, 
1998; Francis & Currie, 2003; Bjorholm et al., 2005), 
while our estimates are based solely on locality data 
from herbarium specimens. These two methods may 
yield different estimates of spatial richness patterns 
and correspondingly different rankings of the impor- 
tance of different hypotheses (Hulbert & White, 2005). 
Because species do not typically occur everywhere 
within the area delimited by range maps (Rondinini et 
al., 2006), richness estimates based on range maps may 
measure richness at a larger grain than estimates based 
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coefficients of the Interim General Model second-generation (IGM2). The dashed lines are isopleths for N Central America 
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on locality data (Hulbert & Jetz, 2007). As such, 
richness estimates based on range maps may tend to 
favor variables representing the SH, related to broad- 
scale species turnover (Hulbert & White, 2005). This 
expectation, however, is opposite to the difference 
between our results and those from previous studies. 
Another important difference is that between 
estimates of plant richness based on floras and 
checklists with no correction for sampling effort 
(e.g., Kreft & Jetz, 2007), and our estimates based 
on rarefaction as an attempt to correct for spatial 
variation in floristic knowledge. Floras and checklists 
differ not only in their geographic extents, but in the 
quality of the data and sampling effort across these 
(Frodin, 2001; al. 2005) Such 
heterogeneity in floristic knowledge can affect 


units Kier et 
estimates of spatial patterns of plant richness and 
estimates of the importance of different hypotheses. 
We attempted to account for this heterogeneity using 
number of specimen records as an estimate of 
sampling effort (Nelson et al., 1990), assuming that 
the number of species found in a sampling unit is a 
function of the number of specimens collected in that 
sampling unit. Nonetheless, it is possible to imagine 
scenarios in which the number of specimens collected 
in a sampling unit is a function of the number of 
species occurring in that sampling unit. We think the 
latter is unlikely because plant species inventories of 
100 X 100 km sampling units across our study area 
are invariably incomplete. Furthermore, the inventory 
for each sampling unit was derived from several 
collecting trips that, together, are bound to obtain 
multiple specimens of species that are common in the 
sampling unit, even if some plant collectors discrim- 
inate against common species. 

Studies also differ in the plant life forms on which 
they focus and, thus, may favor one hypothesis over 
another because various life forms may respond 
differently to different factors (Richerson & Lum, 
1980; Laanisto et al., 2008). For example, tree species 
might be less responsive to spatial heterogeneity than 
1982; Qian & 


Ricklefs, 2004). However, among studies measuring 


nonwoody plant species (Gentry, 


species richness of both woody and nonwoody life 


forms, some conclude that the SE is more important 
than the SH (Francis & Currie, 2003; Moser et al., 
2005; Kreft & Jetz, 2007), while others find the 
opposite (Pausas et al., 2003; this study). It would be 
useful for future studies to specifically address this 
issue by confronting species richness data for different 
life forms against models representing the SE and SH. 

Differences in the central tendency of explanatory 
variables may also explain why our results differ from 
those of other studies in terms of the relative 
importance of the SE and SH. In tropical and 
subtropical regions, where energy input is high (> 
505 mm potential evapotranspiration, Kreft & Jetz, 
2007), plant richness may be largely independent of 
energy and mainly determined by water availability 
(Hawkins et al., 2003; Whittaker et al., 2006; Kreft & 
Jetz, 2007) and, to a lesser extent, spatial heterogeneity 
(Kreft & Jetz, 2007). Despite the fact that potential 
evapotranspiration values were well above 505 mm in 
all of our sampling units (Fig. 3), we obtained 
statistically significant coefficients supporting predic- 
tions derived from the SE about how plant richness is 
related to potential evapotranspiration and minimum 
monthly potential evapotranspiration (Tables 1, 2; 
Fig. 5). Thus, our results suggest that both energy 
and water availability do determine plant richness 
across the study regions but to a lesser extent than 
variables representing the SH. Our sample, however, 
includes few sampling units from some extreme 
environments in the study regions, particularly those 
with the lowest annual precipitation (< 1000 mm 
annual precipitation) and highest water deficit (> 
750 mm; Fig. 3). If the effect of annual precipitation or 
water deficit on plant richness decreases with increas- 
ing water availability (Gentry, 1988; Whittaker et al., 
2003), our results could have underestimated the 
importance of water availability and, therefore, the SE. 

Differences between our results and those of other 
studies in the relative importance of the SE and SH 
may relate to differences in the range of explanatory 
variables. Compared 10 studies spanning global 
extents (Francis & Currie, 2003; Kreft & Jetz, 2007) 
or aiming to develop global models (Field et al., 
2005), the ranges of several variables measuring 


t= 


and the dotted lines for N South America and S Central America, according to the fit of the IGM2 with regional effects to our 
sample (Table 2). —D, E. F. Marginal rates of substitution (MRS) of annual precipitation for logarithm of elevation range (D). 
MRS of annual precipitation for minimum monthly potential evapotranspiration (E), and MRS of the logarithm of elevation 
range for minimum monthly potential evapotranspiration (F). The dashed lines represent 95% confidence intervals of the 
MRSs for N Central America and the dotted lines (and gray shading) for N South America and S Central America. These 
confidence intervals are based on parametric bootstrap samples (Efron & Tibishirani, 1993) of size 1,000,000, assuming a 
multivariate normal distribution of the regression coefficients with the variance-covariance matrix estimated from fitting the 
IGM2 to our sample. The continuous lines represent the MRSs based on the global coefficients of the IGM2 and have no error 
estimates because the variance-covariance matrix for such coefficients was not available. 
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energy and water availability in our sample were small 
(Fig. 3). Therefore, we may have underestimated the 
importance of the SE because the proportion of the 
variation in a response variable (e.g., relative plant 
richness) accounted for by an explanatory variable 
(e.g., mean annual temperature) can be a function of 
the sample range of the latter (Pedhazur, 1997). By 
the same token, in our sample, the range of a variable 
representing the SH, elevation range within sampling 
units, was similar to the sample range of a study of 
global extent (Kreft & Jetz, 2007) and larger than that 
of a study aiming to develop global models (Field et 
al., 2005; Fig. 3). Thus, the latter study may have 
underestimated the importance of the SH. If our 
sample is representative of the variation in energy and 
water availability across our study area (Fig. 3), the 
conclusion of more empirical support for the SH than 
for the SE would still be valid for this area. 

We found significant regional effects among N 
Central America, S Central America, and N South 
America, consistent with studies that have explicitly 
tested for, and commonly found, evidence supporting 
the RE (Schluter & Ricklefs, 1993; Ricklefs, 2004; 
Kreft & Jetz, 2007). Identifying the underlying causes 
of differences in plant richness attributed to regional 
effects can be difficult, because they can actually 
reflect differences in environmental variables not 
accounted for (Schluter & Ricklefs, 1993). However, 
two of the most prominent regional effects that we 
obtained were consistent with previous assessments of 
historical influences on Neotropical plant diversity 
(Gentry, 1982). First, the lower intercept for N Central 
America than for other regions may be due to limited 
of Gondwanan clades that 
compose most Neotropical plant diversity in the 


northward movement 


lowlands (i-e., Amazonian-centered taxa sensu Gentry, 
1982). Second, a steeper slope relating plant richness 
to spatial heterogeneity in N Central America than in 
other regions may be due to limited southward 
movement of Laurasian clades that are more important 
components of montane than lowland Neotropical 
floras (Gentry, 1982). Laurasian clades may then 
contribute more to the increase in richness associated 
with spatial heterogeneity in N Central America than 
in S Central America and N South America. 

Last, we discuss implications of our findings for the 
perception of the spatial pattern of plant richness 
across the western Neotropics. Our analysis predicted 
peaks of relative species richness mostly in topo- 
graphically complex areas (Figs. 2C, 6), consistent 
with a similar analysis for northwest South America 
(Jiménez et al. 2009) and another recent mapping 
effort (Kreft & Jetz, 2007: fig 3b), but contrasting with 
maps showing higher vascular plant richness in 
lowland areas than in the northern Andes (Barthlott 


et al., 2005; Mutke & Barthlott, 2005; Kreft & Jetz, 
2007: fig. 3c, d). The latter maps are similar to the 
pattern predicted by global coefficients for the IGM2 
(Fig. 2D). designed for estimating broad-scale rich- 
ness of woody plants worldwide (Field et al., 2005). 
Below, we suggest a potential explanation for the 
contrast between different plant richness maps of the 
Neotropics, based on differences between the global 
coefficients of the IGM2 (Fig. 2D) and the coefficients 
yielded by fitting the IGM2 to our sample (Fig. 2C). 

The contrast between the map of predicted plant 
richness based on global coefficients of the IGM2 
(Fig. 2D) and that based on the fit of the IGM2 to our 
sample (Fig. 2C) is underlain by a difference in the 
importance of the logarithm of elevation range relative 
1o annual precipitation. This is illustrated by isopleths 
showing combinations of logarithm of elevation range 
and annual precipitation that, according to a given 
model and holding other variables constant, yield a 
constant predicted richness value (Fig. 7A). The slope 
of the isopleths based on our sample are more 
negative than the slope of the isopleth based on the 
global coefficients of the IGM2, suggesting that the 
importance of the logarithm of elevation range relative 
to annual precipitation is higher in the fit of the IGM2 
to our data than in the global coefficients of the IGM2. 
This difference is quantified by the marginal rate of 
substitution (MRS; Caraco, 1979; Brown, 1988) of 
annual precipitation for the logarithm of elevation 
range, which corresponds to the negative of the slope 
of the isopleths in Figure 7A and measures the 
amount of annual precipitation needed to substitute 
a logarithmic unit of elevation range and maintain the 
same species richness in any given sampling unit (see 
appendix S1 note A in Jiménez et al., 2009). The MRS 
of annual precipitation for the logarithm of elevation 
range is higher in the fit of the IGM2 to our data than 
in the global coefficients of the IGM2 (Fig. 7D). 

The foregoing contrast between the maps of 
predicted plant richness (Fig. 2C vs. D) also results 
from differences in the relative importance of 
minimum monthly potential evapotranspiration and 
annual precipitation. This can be seen in isopleths 
showing combinations of values of minimum monthly 
potential evapotranspiration and annual precipitation 
that, holding other variables constant, yield a constant 
predicted richness value (Fig. 7B). These isopleths 
bend because, according to the IGM2, the relationship 
between minimum monthly potential evapotranspira- 
tion and plant richness is quadratic (e.g., Fig. 5A). 
The negative of the slope of these isopleths is the MRS 
of annual precipitation for minimum monthly potential 
evapotranspiration. It measures the amount of annual 
precipitation needed to substitute a small increase in 
minimum monthly potential evapotranspiration and 
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maintain the same species richness in any given 
sampling unit. This MRS decreases as minimum 
monthly potential evapotranspiration increases 
(Fig. 7E) because, as pointed out previously, according 
to the IGM2 the relationship between minimum monthly 
potential evapotranspiration and plant richness is 
quadratic. As minimum monthly potential evapotrans- 
piration increases, the decrease in the MRS of annual 
precipitation for minimum monthly potential evapo- 
transpiration is steeper in the fit of the IGM2 to our data 
than in the global coefficients of the IGM2 (Fig. 7E). 
This finding suggests that the importance of minimum 
monthly potential evapotranspiration relative to annual 
precipitation is higher in the fit of the IGM2 to our 
sample than in the global coefficients of the IGM2. 
Finally, the contrast between the map of predicted 
plant richness derived from global coefficients of the 
IGM2 (Fig. 2D) and that derived from the coefficients 
yielded by fitting the IGM2 to our sample (Fig. 2C) 
does not seem to reflect differences in the importance 
of the logarithm of elevation range relative to 
minimum monthly potential evapotranspiration. This 
notion is suggested by isopleths showing combinations 
of values of logarithm of elevation range and minimum 
monthly potential evapotranspiration that, holding 
other variables constant, yield a constant predicted 
richness (Fig. 7C). These isopleths bend 
because the IGM2 portrays a quadratic relationship 


value 


between minimum monthly potential evapotranspira- 
tion and plant richness. The negative of the slope of 
these isopleths is the MRS of logarithm of elevation 
range for minimum monthly potential evapotranspira- 
tion, and measures the elevation range (in logarithmic 
units) needed to substitute a small increase in 
minimum monthly potential evapotranspiration and 
maintain the same species richness in any given 
sampling unit. Along the axis of minimum monthly 
potential evapotranspiration, the MRS of logarithm of 
elevation range for minimum monthly potential 
evapotranspiration derived from the fit of the IGM2 
1o our data overlaps with that derived from the global 
coefficients of the IGM2 (Fig. 7F). This finding 
indicates that the importance of the logarithm of 
elevation range relative to minimum monthly potential 
evapotranspiration is similar in both cases. 

The comparisons above suggest that differences 
among representations of plant richness across the 
Neotropics may result from differences in estimates of 
the relative importance of three main determinants of 
plant richness. Relative to maps showing higher 
than 
topographically complex areas, our maps of predicted 


vascular plant richness in lowland areas 
richness may assign a larger role to elevation range 
and minimum monthly potential evapotranspiration 


than to annual precipitation. 
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